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Abstract 

Coalescence of low mass compact binaries of neutron stars and black holes are 
primary burst sources for LIGO and VIRGO. Of importance in the early stages 
of observations will be the classification of candidate detections by source-type. 
The diversity in source parameters and serendipity in any new window of ob- 
servations suggest to consider model-independent detection algorithms. Here a 
frequency dynamics algorithm is described which extracts a trajectory in the 
' /(/)"Pl ane from the noisy signal. The algorithm is studied in simulated binary 

coalescence. Robust results are obtained with experimental noise data. Exper- 
j> ■ iments show the method to be superior to matched filtering in the presence of 

^ ■ model imperfections. 

O : 

\o ■ 
o : 

O ■ I. INTRODUCTION 

o ■ 

^: 

This millennium will mark the beginning of gravitational wave astronomy with the broad 
band observatories LIGO |T| and VIRGO 0. This opens up new opportunities to probe 
strongly gravitating processes in the coalescence of binaries and the birth of neutron stars 
or black holes f3|-|5|. The late stages of spiral infall of a neutron star may show interesting 
j> | new physics, for example, in association with tidal break-up around a companion black 
hole, where the onset will provide constraints on the density of neutron star matter. There 
may also be surprises. Advanced LIGO is expected to see events up to a Gpc [HH, and 
might detect correlations or anti-correlations of gravitational radiation with cosmological 
7-ray bursts. The early stages of observations will be focused on verification of potential 
gravitational wave signals and identification of source by type. In view of the diversity 
of source parameters and serendipity in any new window of observations, it is of interest 
to consider model-independent detection algorithms. These algorithms are expected to 
compare well with matched filtering, anticipating the inevitable model imperfections in both 
the wave-forms and the instrumentation noise. 

Early evolution of a binary of two stars of mass Mi and M2 with orbital period P is 
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well- described by quasi- Newtonian evolution in the point-particle limit, whose luminosity 
Cgw i n gravitational wave emissions in the theory of general relativity is given by the Peter 
and Mathews' formula || 

QO 

C GW ~ ~r (ujM) w ^F(e) (1) 
5 

in geometrized units with G — c — 1, where uj = 2tt/ P denotes the orbital angular velocity, 
F(e) accounts for an orbital ellipticity e and M = (MiM 2 ) 3/b / (M x + M 2 ) 1/5 denotes the 
chirp mass. This is accompanied by a gradual rise in the orbital frequency f or b = 1/P: 

L b = ^f*M^flll\ (2) 

The predicted luminosity of 3.2 x 10 33 ergs/s agrees with the observed rate of change 2.40 x 
10~ 12 in the orbital period of about eight hours in the Hulse- Taylor binary system with 
neutron star masses M% = 1.42M Q and M 2 = 1AM & [[/]]. This strong frequency-luminosity 
relationship suggests to plot the evolution of a candidate source as a trajectory in the /(/)- 
plane as the output of a detection algorithm. In doing so, we focus on burst sources whose 
gravitational waves show appreciable dynamics in their frequencies. This, in contrast to 
periodic sources, e.g.: rapidly spinning neutron stars [§,[§, which we regard as a separate 
class of sources in regards to detection algorithms. Binary coalescence of compact objects, 
then, shows an initial branch / oc f 11 ^ 3 indicative of the general relativistic relationship 
Cgw f 10 ^ 3 - These trajectories terminate on the /-axis at the quasi- normal mode frequency 
Jqnr of a final black hole. For the expected low mass black holes M < 1OM in compact 
binaries, /qnr ~ 28OO(1OM /M) is beyond the (advanced) LIGO sensitivity range of about 
(20-1600Hz) 40-800Hz, however. The transition which connects the chirp to the final black 
hole state remains highly uncertain - here the gravitational radiation emission process is 
most nonlinear and potentially most interesting. For black hole-black hole coalescence, the 
transition from chirp to a common horizon envelope state should be smooth, may be very 
luminous Q and more frequent ]10| than neutron star-black hole coalescences [PHT3[| . In 



neutron star-black hole coalescence, on the other hand, an intermediate black hole-torus state 
is expected if the black hole spins rapidly [l4j] . This state is expected to be quiescent in its 
gravitational wave emissions, while luminous in emissions by the black hole in contact with 
the magnetic field of the torus |1J. This predicts an anti-correlation between gravitational 
wave-emissions and 7-ray bursts from this type of catastrophic events |16[ . Alternatively, 
gravitational radiation is expected in accretion induced collapse of a white dwarf into a 
neutron star or black hole due to a bar mode instability [jT7|| , or the collapse of a young 



massive star [18, IS], in which the relationship between / and / is even more uncertain. If 
the collapse stalls with the formation of a neutron star, we might witness a negative branch 
/ < due to spin-down by gravitational wave emissions. It becomes potentially useful, 
therefore, to plot transient gravitational wave emissions as trajectories in the /(/)-plane 
when classifying candidate detections by source- type p0| . As gravitational wave-emissions 



are not expected to be significantly beamed, observations such as these by advanced LIGO 
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can obtain definite evidence of 7-ray burst progenitors and their population statistics which 
are not otherwise available. 

Here we describe an algorithm which enables accurate extraction of the frequency dy- 
namics in the gravitational wave signal. The algorithm is based on counting zero-crossings, 
wherefrom an instantaneous frequency f(t) and frequency rate of change f(t) can be es- 
timated as a function of observer's time t. Here the proposed algorithm is tested using 
simulated binary coalescence with Gaussian noise and noise recorded on the 40m LIGO test 
facility at Caltech. The algorithm is studied on signals of the type 

X(t) = A(t) cos(27r$(t)) + aN(t) (3) 

wth amplitude A(t) = (1 — t) -1 / 4 and phase 

$(t) = c [(l - (1 - tf 8 ) + p (l - (1 - t) 3 / 8 )] , (4) 

where te[0, 1) is the time normalized to the time of coalescence t = 1, c = 1000 and pe[0, 1] 
is a model parameter: p = in the Newtonian approximation and p = 0.6038 in the 



1PN approximation when Mi = M2 pl |. The noise process N(t) is normalized such that a 
controls the signal-to- noise ratio. All simulations are for data strings of a total of N = 10, 000 
samples at t n = nAt (n = 1, 2, • • • N). The instantaneous frequency fit) is given by d<&{t)/dt 
and d 2 $(t)/dt 2 denotes its rate of change fit) with respect to time. The dimensionless rate 
of change of the period, f~ 2 f, expresses inversely the number of periods at a frequency 
/, which is of the order of a few hundred in the experiments reported below (final time 
t ~ 0.97). 



II. DESCRIPTION OF THE ALGORITHM 

Consider the instantaneous frequency f(t) as a finite time-series for t e [a,b]. We begin 
with a partition of [a, b] into N adjacent subwindows (bins) [sj_i,Sj] with centers at U = 
(si_i + Si)/2. Linear approximations f(t) ~ a« + (3$ (sj_i < t < Si) are then obtained by 
estimates of the coefficients a» and $ from the noisy data X(t) in each subwindow i. The 
number of periods in the z-th subwindow, for example, is f** f(t)dt = a^8 + (3iU5 + 0(5 3 ), 
where S — — s,_i. 

Let Zs(si-i,Si) denote the number of zero-crossings of the true signal S(t) in the sub- 
window [sj_i,Sj]. Then Zs{si-i, Si) differs from 2 f**_ f(t)dt by at most one. In view of 
f(t) = d<&(t)/dt, our algorithm hinges on the approximation 

Z s (s i - 1 ,s i )~26a i + 26l3iti. (5) 

The number of zero-crossings in the noisy data (X (uo) , X (u\) , ■ ■■ ,X(uk)) in will 
be denoted by s i), where the Uj denote the discrete sampling times provided by 

the AD-converter (16 bit resolution, 10kHz). The X(uj) have expectation values S(uj) and 
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variance a 2 . The lOkHz-sampling frequency is assumed to be sufficiently high such that the 
true signal agrees with the zero-crossings of its linear interpolant based on the S(uj). In 
terms of the indicator function, /, defined by I(P) = 1 if P is true and I(P) = if P is 
false, we have 

Z s ( Sl . u Si ) = E^/[5(^_!)5(^) < 0], Z x ( Sl . u 8i ) = S^[X(^_ 1 )X(n,) < 0]. (6) 

The probability pj that the observed signal X(uj) and the its expectation (the true signal) 
S(uj) have the same sign can be easily computed. Suppose that S(uj) > 0. Then, we have 

P[X( Uj ) > 0] = P[(X( Uj ) - S{u 3 ))/a > -S( Uj )/a] = 1 - F[-S(u 3 )/a] = F[S( Uj )/a], (7) 

where F is the probability distribution of the noise process N(t), which is assumed to be 
symmetric about zero. Similarly when S(uj) < 0, we find P[X(v,j) < 0] = F[—S(uj)/a]. It 
follows that pj = F[\S(uj)\/a], whereby 

P[I[X( Uj -i)X( Uj ) < 0] ^ I[S( Uj ^)S( Uj ) < 0)]] =p j -i{l-p j )+p j (l-p j -i). (8) 

Consider, then, Ri = mmj\S(uj)\/a, as a measure for the signal-to-noise ratio in the sub- 
window [sj_i,Sj]. As Ri — > oo, miiijPj — > 1. We then have for the expectation E and the 
variance V the limits 

lim J5[Zx(si-i,Si)] = Z s [si-i,Si), V[Z x {8 i -x,8 i )\ = Q. (9) 

Ri— »oo 

Hence, for large signal to noise ratios, we recover Zx(si-i, s») ~ 25cti + 28(3^. 

Estimates «j for f(ti) and /3j for /(tj) can be calculated by linear regression. Define (as 
is usual in such regressions, see, e.g. ||22|| ) the L 2 error 

Qfa, A) = Sf 8i ) - (25ai + 25f3 l t l )) 2 . (10) 

Let &i, (3i minimize Q(a,f3). The regression estimates of the instantaneous frequency f(ti) 
and its rate of change /(tj) are and /3j respectively. As these are obtained simultaneously, 
the / are on equal par with /, rather than being their derivatives. Standard regression 
theory allows us to compute error boxes in these estimates. Note that the method amounts 
to an approximation of the time-series f(t) by a piecewise linear, and generally discontinuous 
function. 

The method may be optimized by modifying the present L 2 norm in the error function 
to curve fitting, e.g.: different weights on / and /. Initial experiments indicate that the 
method can be implemented with adaptive window sizes, which may improve its efficiency 
when the dynamic range in the frequency f(t) is large. 

Because only the pattern of zero-crossings is measured, the algorithm targets precisely the 
underlying gravitational wave frequency-luminosity relationship of the source, as examplified 
in binary coalescence. The algorithm suppresses the amplitude of the signal, and hence it 
is robust to size fluctuations in the signal. Likewise, reference to the underlying cosine 
functions is weak. In fact, any 27r-function with two roots in each period can replace the 
cosine function in (1). 
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III. NOISE DEPENDENCE 



The frequency dynamics algorithm has been tested in simulated binary coalescence using 
the Newtonian model of spiral in (p = 0) with Gaussian noise and the instrumental noise 
from the 40m-LIGO test facility at Caltech. The chirp extends over t = [0,0.97), wherein 
the number of periods at a given frequency f(t) decreases from above 400 to below 200. 
Figure 1 shows the results for both cases. 

The algorithm obtains higher accuracy with the 40m-LIGO noise than with Gaussian 
noise when the variance is high. We believe that the reason for this lies in the high-degree 
of dependence in the noise from the LIGO prototype. The first-order autocorrelations from 
40m-LIGO noise were of the order of 0.8. There were, on average, 15 zero-crossings per 100 
recordings of the pure noise process. This is much smaller than the 50 or so zero crossings 
obtained for white noise of the same length. Thus the LIGO noise contributes far fewer 
"false crossing" due to noise as opposed to signal. This will be the case for all noise with 
low-frequency components in the spectrum. In such cases, our algorithms performs even 
better than it would with white noise. To make the last point more precise, consider a 
Gaussian process with a first order autocorrelation of p. If we observe a sequence of length 
N from such a process the expected number of zero-crossings is 



For p of 0.8, as in the LIGO case, we would expect approximately 20 zero crossings per 100 
recordings. Of course, the LIGO noise from the Caltech prototype is not pure Gaussian, but 
the principle is the same. 



Optimal application of matched filtering requires accurate models of the system to be 
detected and the noise of the detector. Both aspects of prior knowledge are critical and 
must be described in a low number of parameters, lest matched filtering becomes unstable. 
Detection of gravitational waves from forementioned candidate sources poses a number of 
uncertainties which may enter as model imperfections (e.g., size of neutron stars [0]). These 
imperfections can give rise to woefully wrong results when the true signal is not contained 
in the model for any choice of its parameters. It becomes of interest, therefore, to compare 
the model-independence aspect of the frequency dynamics algorithm with the sensitivity 
of matched filtering to these imperfections. The comparison would be complete with an 
additional on noise type and intensity which, however, falls outside the scope of this work. 
The comparison is made using the post-Newtonian model in the 1PN approximation as the 
true signal (p = 0.6038), while applying matched filtering assuming the Newtonian model 
given by p = 0. This mismatch of model to true signal is used here to simulate other model 
imperfections, such as due to the finite size of the objects. 




IV. COMPARISON WITH MATCHED FILTERING 
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Thus, the matched filtering model is perfect in case of p = 0, and becomes imperfect 
as p 7^ 0. It is expected that matched filtering will outperform the frequency dynamics 
algorithm for p — 0, and vice-versa for certain larger values of p > 0. The comparison is 
made quantitative using a mean absolute relative error (MARE), expressed as a percentage, 
in the estimated frequencies by comparison with the known signal input: 

100 \f(ti) - f(U)\ , x 

Using a small Monte Carlo study, the comparison is made at three levels of additive Gaussian 
white noise. Figure 2 shows the results of this comparison study, where each reported MARE 
is the average of 200 simulations. 



V. CONCLUSIONS 

The trajectory in the /(/)-plane of a burst of gravitational radiation is expected to 
provide a signature of the source. This is anticipated from the strong - possibly characteristic 
- coupling between the frequency and the luminosity of its gravitational wave emissions. 
Binary coalescence, for example, shows the characteristic chirp / oc f a+1 / 3 associated with a 
luminosity Caw oc f a , where a = 10/3 in general relativity in the point-particle limit. Other 
burst sources of gravitational radiation may have different frequency dynamics (e.g.: different 
en's). We have described a nonparametric (model-free) frequency dynamics algorithm to 
obtain an accurate plot of the trajectory in the /(/)— plane. It is based on zero-crossings 
in the gravitational wave signal, to extract these trajectories in the presence of noise. The 
method is robust against non-Gaussian additive noise and amplitude variations. The method 
shows good results at appreciable chirp rates and different types of noise. Even when 
the signal-to-noise level is low (down to about 2 in our experiments), we can still recover 
the distinctive curve of the frequency dynamics when / is high. The estimation appears 
to be much more accurate for LIGO noise than for Gaussian noise when the variance is 
high. We believe that the reason for this lies in the high-degree of dependence in the 
noise from the LIGO prototype, which suppresses "false zero-crossings" generated by the 
noise. This will be the case for all noise with low-frequency components in the spectrum. 
In such cases, our algorithm performs even better than it would with white noise. The 
frequency dynamics algorithm has been compared with matched filtering, and shown to be 
superior in the presence of model imperfections. This is particularly relevant in searches for 
unanticipated burst sources. 

The present experiments indicate that the algorithm should perform well in extracting 
the frequency dynamics in the gravitational waves. Future experiments are desirable at very 
low signal-to-noise levels, to study the detection limit for transient sources with the present 
algorithm. 
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Figure Captions 



Figure 1. The /(/)-diagram of simulated binary coalescence in the Newtonian approxima- 
tion, extracted by the frequency dynamics algorithm in the presence of noise at a signal-to- 
noise ratio of about 4 (s = a = 0.25). The horizontal axis shows the frequency / [Hz] and 
the vertical axis its time rate of change df [Hz/s]. The top window shows the results with 
Gaussian noise, and the lower window shows the results with the instrumental noise from 
the 40m-LIGO test facility at Caltech. The crosses in the right windows denote horizontal 
and vertical errors obtained from the linear regression estimates. The results indicate that 
the method improves when the noise is correlated. At higher frequencies, the results are 
similar due to a larger instantaneous signal-to-noise ratio. 

Figure 2. Experiments on model-independence in terms of the normalized error MARE on 
the frequency, obtained by the frequency dynamics algorithms (lines) and matched filtering 
(dashed) at three different variances a 2 = 0.01,0.15,0.25. The comparisons are made for 
three different signals, parametrized by the value of p (horizontal axis). Matched filtering 
is performed with model assumption p — 0, hence p ^ in the signal simulates a model 
imperfection. The results show that the frequency dynamics algorithm outperforms matched 
filtering in the presence of such appreciable model imperfections, with MARE approaching 
a constant. 
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Figure 2. 
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